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Three-dimensional Ising model confined in low-porosity aerogels: 
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The influence of correlated impurities on the critical behaviour of the 3D Ising model is studied 
using Monte Carlo simulations. Spins are confined into the pores of simulated aerogels (diffusion 
limited cluster-cluster aggregation) in order to study the effect of quenched disorder on the critical 
behaviour of this magnetic system. Finite size scaling is used to estimate critical couplings and 
exponents. Long-range correlated disorder does not affect critical behavior. Asymptotic exponents 
differ from those of the pure 3D Ising model (3DIS), but it is impossible, with our precision, to 
distinguish them from the randomly diluted Ising model (RDIS). 
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I. INTRODUCTION 

The influence of quenched disorder on phase transi- 
tions has been studied since long time ago now. In 1974, 
Harris^ established this famous criterion: uncorrelated 
disorder is not relevant, for a second order phase tran- 
sition, if the specific heat exponent is negative (a < 0). 
The criterion was generalized by Weinrib and Halperin 
(WH)^ for any disorder distribution whose correlation 
function exhibits a power-law tail, i. e., g{r) ~ r"" as 
r — > oo. Disorder is shown to be relevant in these cases: 



di^ - 2 < if a > d, 
ai^ - 2 < if a < d, 



(1) 
(2) 



being d the dimension, and v the correlation length ex- 
ponent of the pure system. After Josephson hyperscaling 
(2 — di/ = a) Harris criterion is recovered in the short- 
range correlated (SRC) regime (0. In contrast, the long- 
range correlated (LRC) regime extends the criterion to 
systems satisfying condition Q, even if a < 0. This 
generalization explains why critical exponents for the su- 
perfluid (SF) transition of '*He change when the fluid is 
conflned in very light aerogels'^''* and do not when con- 
fined in, for example, porous gold^. Aerogels are fractal 
for several length scales-:, while porous gold has exponen- 
tially decaying correlations beyond the size of a typical 
poreSi. Nevertheless, authoraM argued that the critical 
behaviour of SF ^He in aerogels yet poses intriguing ques- 
tions to be solved. 

For instance, light aerogels are fractal for several length 
scales, up to a certain value A that depends on aero- 
gel density. Beyond this length, the disordered structure 
becomes homogeneous, entering an uncorrelated regime. 
After Harris criterion, as the SF correlation length ^ 
gets larger than A, disorder should become irrelevant, 
because a ~ —0.011 is negative for this system. Yoon 
et 

t : , 
Tc as close as t ~ 10^^, no crossover to bulk exponents 



. estimated that this crossover should appear at 
T - Tc\/Tc ~ 10""* but, although they approached 



was observed. A different universality class was evident 
for the SF transition of *He, when confined in aero- 
gels. An explanation to these changes was given using 
Monte Carlo (MC) simulations of the three-dimensional 
XY (3DXY) model, confined in aerogel-like structures^. 
The SF transition belongs to the 3DXY universality class, 
and correlated disorder could be relevant provided that 
the WH condition (0) at r ^ cx) is fulfilled. Vasquez et 
al7'^ showed that changes occur because of hidden LRC, 
inherent to the process of aerogel formation. Using sim- 
ulated aerogels, made by diffusion limited cluster-cluster 
aggregation (DLCA)-'-- , authors showed that different 
LRC subsets are physically well defined within the whole 
aerogel structure. Specifically, gelling clusters (to be de- 
fined later in this paper) are shown to be the relevant 
structures defining the critical behaviour of the 3DXY 
model in DLCA aerogels^. 

In this paper, we study the three-dimensional Ising 
model (3DIS), in presence of such aerogel-like structures. 
The pure 3DIS has a positive specific heat exponent (a ~ 
0.11), so any type of disorder, correlated or not, will be 
relevant. If Ising spins are collocated in the pores of 
aerogels, criticality will be affected by LRC as well as by 
SRC disorder. Our main purpose is to elucidate which 
among these effects dominates the critical behaviour of 
the 3DIS model in this case. Along this paper, we report 
the results of extensive MC simulations of the 3DIS in 
the pores of DLCA aerogels at fixed porosity ip = 80%, 
in order to clarify this point. 

The rest of this paper is organized as follows: Section 
^ is a brief review about diluted Ising systems studied 
in the past. Section IlIII is due to explain the model first, 
then the simulation procedure in detail, with a pream- 
ble on self-averaging, in order to validate our procedure; 
two methods to obtain accurate values of the critical cou- 
pling are presented at the end. Thermal and magnetic 
effective exponents are presented, and their asymptotic 
behaviours are discussed in Section llVI Finally, section 
fVlis due to present some concluding remarks. 



II. ANTECEDENTS ON THE DILUTED 3DIS 

Many experimental, theoretical, and computational 
works have been done to date, in order to study the 
critical behaviour of the 3DIS model in presence of 
quenched disorder. Most of numerical and theoretical 
works address the random-diluted 3DIS (RDIS), i. e., 
the Ising model in presence non- correlated distributions 
of impuritiesiiiiSiiii^ii^. 

Appart, concerning disordered LRC structures, Balles- 
teros and ParisiiS- simulated the 3DIS, with dislocations 
represented by lines of impurities generated at random. 
Correlations for this type of disorder decay with an ex- 
ponent a — 2. They obtain a correlation length exponent 
i^LRc ~ 1 for the impure system, thus confirming the re- 



sult of WH, that this exponent should be 



2/a. 



Marques et alm^ also simulated a diluted 3DIS, but with 
spins located on LRC sites. These sites were provided by 
a previous simulation of the pure 3DIS model; taking all 
sites from the cluster of dominating spin orientation at 
Tc, these are then occupied by the interacting Ising spins 
to simulate. Clusters for this thermally diluted Ising sys- 
tem have anomalous dimension rj^^^^ « 0.03, which gives 
a = 2 — 77pujo ~ 1.97. They obtain an exponent i^lrc ~ 1, 
also in agreement with the WH expression. In both cases, 
LRC disorder is relevant for criticality. Nevertheless, this 
particular result from WH has proven recently not to be 
correct at more accurate approximations. Using two- loop 
expansions, Prudnikov et alm^, showed that the exponent 
t'LRc depends on both, the internal dimension of the or- 
der parameter m and the exponent a, not the case in 
WH's conjecture, independent of m. For both systems, 
Prudnikovs' prediction yields i'lrc ~ 0.72 ^ 1. 

Experiments about the critical point of the liquid- 
vapour (LV) transition of ^He and N2, confined in 95% 
porous aerogels^^'^", concern directly the problem we are 
addressing in this paper. Bulk '^He near its LV critical 
point belongs to the 3DIS universality class, and aerogel- 
like disorder has proven to contain both, LRC and SRC 
disordered structures^iSi. Wong et ali^^ report, for the 
order parameter, exponents (3 = 0.28(5) and 0.35(5), re- 
spectively, consistent with that calculated for the pure 
3DIS by Guida and Zinn-Justin2i (3 = 0.326(3). How- 
ever, in the same experiments, the specific heat curves 
present finite peaks at Tc, characteristic of a negative 
exponent a, definitely different from the corresponding 
Qipure ~ 0.11 for the pure 3DIS. Actually, within error 
bars, results for /3 are also consistent with the correspond- 
ing RDIS valueS^, /3 = 0.355(5). After these experimen- 
tal results one may take non-correlated, SRC instead of 
LRC disorder, within aerogels, to be the relevant one for 
the critical behaviour of the 3DIS in aerogel pores^. 

Renormalization group (RG) calculations for the 3DIS 
with weak amounts of disorder show that a new uni- 
versality class appears, different from that of the pure 
3Pjgi2^i3,24^25,26,27,28^ g^^^ Consistent with Harris crite- 
rion. Since D. P. LandauSS, using MC simulations, con- 
cluded that the exponents for the 3DIS with random im- 



purities differ from those of the pure system, different 
works stated that exponents depend on the concentra- 
tion of impurities. Until 1990, when Heueri^Si^ began to 
clarify that differences with RG calculations are due to 
the fact that exponents obtained from simulations were 
basically effective ones, and not the asymptotic ones. 

Ballesteros et alm^, using a p-reweighting method in 
MC simulations, found the exponents for the RDIS uni- 
versality class to be independent from the concentration 
of impurities p. This was confirmed for the random bond 
Ising model in d = 3 by Berchc et al.'^^, looking at finite 
size scaling for the critical temperature. All those cal- 
culations motivated further MC and RG studies about 
the crossover between the effective and really asymptotic 
critical behaviourii^iSSiS. 

Definitive evidences of a new universality class were ob- 
tained using neutron scattering in the antiferromagnets 
Mni^2;Zna;F2^ and Fe^Zni-ajFj^. For the first system, 
exponents v = 0.70(2) and 7 = 1.37(4) are obtained, 
while v = 0.69(2) and 7 = 1.31(3), are the results for 
the second one. These results clearly differ from bulk ex- 
ponents, see for instance those calculated by Guida and 
Zinn-Justin22: j/ = 0.6304(13), 7 = 1.239(5). 

The critical behaviour of magnetic systems confined 
in aerogel-like structures, may be subject of competing 
LRC and SRC influences. The 3DXY model in pores 
of DLCA aerogels, for instance, presents new exponents 
due only to the presence of gelling clusters, which are 
LRC, while the SRC components are irrelevant to the 
transition^. However, for the Ising model under the same 
kind of confinement, two different effects may be present. 
Simulations under strictly LRC types of disorderiiSiii, 
give exponents consistent with the result of Weinrib and 
Halperin. On the other hand, experimental results about 
the critical point of LV transitions in aerogelsiSiSS point 
to the relevance of the uncorrelated part of disorder. 



III. THE MODEL AND SIMULATION 
PROCEDURE 



The 3DIS in the presence of impurities, on a simple cu- 
bic lattice with nearest-neighbor interaction, is described 
by the Hamiltonian: 



/3W = -T— ^e.ejs.sj-. 



(3) 



where s = ±1 are the spin variables, fcg is the Boltzmann 
constant, and J is the coupling. In what follows, kg = 1 
and T — 1. The sets {ei} represent quenched variables 
chosen to be if the site is an impurity and 1 if the site 
is occupied by a spin. These sets of impurities are taken 
randomly (RDIS) or from DLCA aerogels (AEIS). 



A. Disorder generation and MC simulations 

At the beginning, sites are occupied by a uniform ran- 
dom distribution of N particles, so their volume fraction 
is c = N/L^. To simulate the RDIS, this initial distri- 
bution of disorder is held through the rest of the MC 
simulation. 

Instead, for AEIS simulations to take place, aerogels 
are generated through the on-lattice DLCA algorithnAiS 
with periodic boundary conditions (PBC). Monomers 
and clusters diffuse randomly with diffusivity constants 
D which depend on their mass n through D ~ n~^''^f . 
The fractal dimension df has been taken equal to its value 
in three dimensions'^'^, df « 1.8. They stick irreversibly 
when they come in contact, and then the process fol- 
lows up until a single cluster is obtained. This model is 
known to reproduce well the geometric features of real 
aerogels^k. 

At a given stage in this DLCA process, the first ag- 
gregate to reach opposite sides of the simulation box in 
any direction is called the gelling cluster (GC). It has 
been shown£*Si that the correlation function for this GC 
is algebraic up to a cutoff, which diverges as L — > cx). 
In other words, these objects are fractal (LRC). Right 
after the GC is built, many other smaller clusters (is- 
lands) continue to diffuse, and finally attach themselves 
to the GC at random sites. The resulting DLCA cluster 
(GC with islands) becomes homogeneous at a very small 
cutoff, in spite of the existence of a physically well de- 
fined fractal structure, the GC. This cutoff increases as 
the concentration decreases, a feature already observed 
for real silica aerogels^. At the volume fraction employed 
in the present work, c — 0.2, the cutoff is so small (a 
few lattice constants) that DLCA clusters must be con- 
sidered as non-fractals. Thus, the presence of islands, 
which represent the SRC subset within the whole DLCA 
cluster, actually hide the LRC behaviour of gelling clus- 
ters. It is in this sense that aerogel-like structures must 
be considered as a mixture of LRC and SRC disorder 
distributions. 

Disordered samples are generated by the procedures 
described above. A MC simulation is performed for 
each sample of interacting 3DIS spins, placed at empty 
sites left by impurities. Physical observables, denoted 
by caligraphs O, are calculated at each indepedent MC 
step, and then corresponding ensemble (thermal) av- 
erages (O) are taken over the MC time-series. Wolff 
algorithm^ is used to update spins. In disordered sys- 
tems, this algorithm tends to prevent some regions from 
being visited by growing Wolff clusters. If the concen- 
tration of impurities is small, this problem can be solved 
by adding some Metropolis updates along the simula- 
tion process.2i. We chose this method and include some 
Metropolis sweeps to shake all spins, after a fixed num- 
ber of Wolff steps. An independent step is taken after 
one correlation time r, which has been estimated from 
preliminary simulations. After enough steps for ther- 
malization, a fixed number (Nt = 1000) of indepen- 



dent MC steps are performed to calculate thermal av- 
erages. Equivalent simulations take place for Ng = 2000 
different samples and, finally, averages over disorder are 
taken O = [{O)] (denoted by square brackets). System 
sizes are L = 8, 12, 16, 24, 32, 48, 64, 96 for the RDIS, and 
L = 8, 12, 16, 24, 28, 32, 40, 48, 56, 64, 80, 96 for the AEIS. 



Measured observables. 



The magnetization (order parameter) is calculated by 



M 



T ■= 

-E 



^i-^ii 



(4) 



where N = cL^ is the total number of spins. Ther- 
mal averages (M) are taken, and averages over disorder 
M{J) = [(Ai)] (J), are then calculated after the former 
have been extrapolated by reweighting'^^. The procedure 
is described below in detail. In terms of the magnetiza- 
tion, we define the susceptibility as: 



X = JL^ [{M^) ~ {Mf] . 
The energy is correspondingly defined by 

£ = -J^eiEj-SiSj, E = [{£)] 



(5) 



(6) 



and then the specific heat is obtained from fluctuations 
of the energy: 



c. = L-' [{£') {£f] 



(7) 



Logarithmic derivatives of rt*'* moments Af" of the 
magnetization [n — 1,2,4), respect to the coupling, are 
calculated through the energy-magnetization covariance: 



91n(X") 



dJ 



(7W"g)-(g)(7W^ 
(TW") 



(8) 



C. Simulation temperatures and rew^eighting. 

Vasquez et al7 found the phase diagram for the 3DXY 
model in the pores of DLCA aerogels. They obtained 
Tc{c)/Tc{0) = Jc{0)/Jc{c) as a function of the con- 
centration of impurities c, being Jc{c) the 3DXY crit- 
ical coupling at volume fraction c of the aerogel. The 
shape of this phase diagram comes basically from the 
porous structure of disorder, specially at low concentra- 
tions. Using this information and the critical coupling 
Jc{0) = 0.2216595(26) for the pure 3DISH, Vasquea^ 
made a rough estimate for the critical coupling for the 
3DIS in the pores of DLCA aerogels at c = 0.2. Mak- 
ing simulations at this rough estimate, and using lat- 
tice sizes L = 10 — 80 and finite size scaling, the value 
Jc(0.2) = 0.25855(3) is obtained for the critical coupling. 



Although for those simulations the number of disorder 
realizations is low (iVg — 30), they obtain critical expo- 
nents close enough to those reported for the RDIS. 

In the present work, all simulations were done at 
c = 0.2, using simulation temperatures Jc = 0.285745 
for the RDIS (following Calabrese et aZ.-^), and at the 
above estimate J^ — 0.25855 for the AEIS. Physical 
quantities at J ~ J^,, are obtained by the reweighting 
method introduced by Ferrenberg and Swendsen^^. This 
procedure was used for each disorder realization at each 
system size L. Each thermodynamic quantity was then 
averaged over disorder for each J within the extrapola- 
tion interval. Finally, maxima of [x(J, L)], \c(^J^Ly\ and 
[91n(A^")/9J(J, i)] were obtained from averaged curves, 
with their corresponding pseudocritical couplings J*{L). 



D. Disorder sampling and self-averaging. 

This part is dedicated to determine a suitable number 
Ns of disorder realizations to obtain critical exponents 
for the AEIS at enough accuracy. A complete study on 
probability distributions for different thermodynamic 
quantities is performed, and results compared with the 
corresponding ones for the RDIS to decide Ns- As an 
example, in Fig. Q] we depict critical susceptibility points, 
obtained for corresponding disorder samples in the RDIS 
(left) and AEIS (right) cases (both at c = 0.2). Points 
come from simulations at the couplings Jc estimated 
above, using the largest lattice sizes (L — 96). This 
distribution looks sharper in the AEIS case, and notably 
more symetric respect to the average than in the RDIS 
case. Top of Fig. [5] shows the probability distribution 
for the susceptibility in both cases. 
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FIG. 1; (Color online) Distribution of the susceptibility for 
the RDIS (left) and AEIS (right) types of disorder, at a con- 
centration of impurities c = 0.2 and lattice size L = 96. Sim- 
ulations were performed at Jc — 0.285745 for the RDIS—, 
and at Jc = 0.25855 for the AEIS. Running averages over the 
samples [xi] are shown by black thick solid lines. 
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FIG. 2: (Color online) Probability distribution of the nor- 
malized susceptibility x — Xi/ [xl- (Top) Distribution for the 
AEIS looks sharper and more symetrical than that for the 
RDIS (data taken from Fig. . (Bot.) Distributions at differ- 
ent lattice sizes for simulations of the AEIS. The distribution 
width appears independent on the lattice size. 

As remarked, the distribution for AEIS (filled circles) 
is sharper and more symetric than the distribution for 
RDIS (empty squares). Note that, in terms of the nor- 
malized variable x — Xi/ [x]j the maximum of the dis- 
tribution for AEIS is closer to the average x — 1 than 
for RDIS. Bottom of Fig. |21 shows the probability distri- 
bution of susceptibility for different lattice sizes. Inde- 
pendence of the width of distributions from L is clear, 
typical for systems lacking of self-averaging, which is the 
expected behaviour of thermodynamic quantities for any 
disordered system at criticality^. 

Self-averaging can be quantitatively checked by the 
normalized squared width-*^ Ra'- 



Ra{L) 



[A^L)] [A{L)f 
[A{L)f 



(9) 



being A any given thermodynamic quantity. In this 
paper, Ra{L) was estimated for the RDIS to com- 
pare with previously reported values. We obtain, as 
L —> cx), Rm — > 0.054 for the magnetization, and 
R^ -^ 0.016 for the susceptibility, both in agreement 
with previous results^^. The ratio here obtained 
Rm/Rx — 3-4 disagrees with RG predictions: Aharony 
and Harris-'' obtained, using e = 4 — d expansions, that 
the leading term is Rm/Rx = 1/4- The discrepancy 
may come from higher order terms in the expansion, 
and not from the definition of the susceptibility as 



was suggested by Berche et alm^. Note also that, in 
the present work, the definition for the susceptibihty, 
X = JL^[{M^)-{M)^], differs from that used by 
Wiseman and Domany"'^, y = JL^ \(A4'^)] . 



as the values of the normalized squared widths Ra are 
substantially lower than those obtained for the RDIS. 



E. Critical couplings 
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FIG. 3: (Color online) Normalized squared widths Ra ver- 
sus inverse system size L~^ for the AEIS at criticality. [Top) 
Asymptotic non-zero values for Ra as L — » oo, for the mag- 
netization (M), the susceptibility (x), and the specific heat 
(ch), evidence the lack of self- averaging in these quantities. 
{Bottom) The power-law behaviour Re ~ L~^ , with a fit- 
ting exponent x ~ 2.58 < d, indicates that the energy {E) is 
weakly self-averaged. 

Results for Ra{L), plotted versus L^^, in the AEIS 
case, are shown on Fig. 13 being A the order parameter 
M, the susceptibility x, the specific heat Ch, and the en- 
ergy E, respectively. As observed for M, x and Ch (top), 
Ra tends to non-zero values as L — > cxd, though asymp- 
totic limits for Ra are smaller for the AEIS than for the 
RDIS, as expected {Rm -^ 0.020 and R^ -^ 0.0028, as 
L — > oo). The power-law behaviour Re ~ L~^ has been 
depicted for the energy (bottom) , being the fitting expo- 
nent X ~ 2.58. Thus, the energy is weakly self-averaged^ 
{x < d). The same type of behaviour is obtained, in 
this work, for the energy in the RDIS case, in agree- 
ment with previously reported results^^. According to 
these analyses, the number Ns of disorder realizations, 
suitable to estimate critical exponents, is larger for the 
RDIS than for the AEIS. In next section, we report some 
values for effective critical exponents for both models. 
Our results for the RDIS agree well with previously pub- 
lished results^^i, thus, the same number of realizations for 
the AEIS will be enough to estimate critical exponents. 



We use two methods to estimate the critical coupling 
out of our present simulations. Binder magnetization 
fourth cumulant 



C/4 = l- 






(10) 



is universal, i. e., independent on the system size^ at 
the critical point. Thus, the critical coupling Jc can 
be obtained with high accuracy at the point where C/4- 
J plots coincide for all system sizes L. Fig. ^ shows 
these plots for the largest system sizes, for the RDIS 
(top) and the AEIS (bottom). In both cases, the cou- 
pling used in simulations has been marked by a verti- 
cal dashed line. For the RDIS case, taking the curves 
for the two largest sizes (L = 64,96), we obtain an esti- 
mate (circle) Jf ^^"^ = 0.2857471(11) for the critical cou- 
pling, which agrees well with that of Calabrese et aLl^, 
jRDis ^ 0.2857447(24). In the AEIS case, we depict 
the same plots for lattice sizes L = 56 — 96. The inter- 
section for L — 80, 96 (circle) gives the critical couphng 
jAEis ^ 0.258575(10), close to the value used in simu- 
lations, J^^^^ = 0.25855 (AJc/Jc ~ 10-''). 

On the other hand, following the FSS theorji^, devia- 
tions of pseudocritical couplings J* (L) from the critical 
coupling Jc scale as: 



j;(L)-j, ^L-l/^ 



(11) 



where Jc{L) is defined as the positions of maxima for a 
given critical quantity, being Jc = Jc{L ^ 00). For in- 
stance, values J*(i) for the susceptibility Q and loga- 
rithmic derivatives © , obtained from reweighted curves, 
have been depicted in Fig.Oas functions of L"^/"^, where 
the rough value l/v « 1.4 has been estimated through 
non-linear fits of points corresponding to L = 56 — 96. 
As expected, the linear behaviour ((TT|l is observed, and 
the lines cross the L — > cx) axis at an average point 
Jc ~ 0.258570(13), quite close to the value estimated 
above, using the intersection of Binder fourth cumulants. 



IV. EFFECTIVE EXPONENTS 
A. Correlation length exponent 

Finite size scaling (FSS)-*^ has been used to estimate 
effective critical exponents for the RDIS and the AEIS. 
This method allows us estimate critical exponents P/v, 
^ /v, otjv-, and the correlation length inverse exponent 
Xjv. The latter has already been roughly estimated 
above, using the scaling law for the position of the 
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FIG. 4: Binder magnetization fourth cumulant as a function 
of J for the RDIS (top) and the AEIS {bottom), using the 
largest L. Dashed vertical lines indicate couplings Jc used in 
simulations, empty circles are intersections between curves at 
the two largest sizes. In the AEIS case, Jc from this intersec- 
tion and the previous estimate Jc differ by AJc/Jc ~ 10~*. 
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FIG. 5: (Color online) Position of the maxima J* for the 
susceptibility and logarithmic derivatives of M" (n = 1, 2, 4), 
plotted versus L~^' "^ using 1/v — 1.4. Linear fits average at 
J = 0.258570(13) in the thermodynamic limit (L -^ oo). 



maxima of logarithmic derivatives of the magnetization 
moments'^'* {A4'"-) {n = 1,2,4), and the susceptibility. 
More accurate estimations are made directly taking aver- 
ages over disorder on quantities obtained at Jc = 0.25855, 



used in our extensive simulations, which is close to the 
previously estimated Jc'- 



dln{M") 



dJ 



il/^ 



(12) 



j=Ja 



We look first to this exponent in order to determine the 
effect of disorder on the critical behaviour of the 3DIS. 
Previous works'^'^ report that exponents fi/iy and j/v for 
the RDIS are almost the same as those for the pure 3DIS, 
and as shown later, this is the case for the AEIS. 

Logarithmic derivatives of moments n = 1,2,4 of the 
magnetization arc plotted, versus system size L, in Fig. El 
for the AEIS case. Points were obtained from averages 
over disoder at the simulation coupHng Jc — 0.25855 
which is quite close to J^ as estimated above. Dashed 
lines are power-law fits to the equation (|12|l using the 
four largest system sizes (L = 56 — 96), and give a FSS 
exponent 1/v = 1.501 ± 0.007 for the AEIS. We recaU 
that exponents determined by this method are effective 
exponents, and only their asymptotic behaviour would 
give a hint to what the universal critical exponent tend 
to. This study is addressed in next paragraphs. 
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FIG. 6: Logarithmic derivatives [dln{M")/dJ] (n = 1,2,4) 
versus L. Points are averages over disorder taken from our 
simulations at the estimate Jc = 0.25855 for the AEIS. Power- 
law fits (dashed lines) out of the four largest lattice sizes L — 
56 - 96 give a FSS exponent l/u = 1.501 ± 0.007. 

Effective exponents (l/i^)eff have been depicted in 
Fig. [3 as calculated from FSS of logarithmic derivatives 
in both cases, RDIS and AEIS. As in Fig. El values were 
obtained from averages over disorder at the simulated 
couplings Jc ~ Jc, being Jc = 0.285745 for the former, 
and Jc — 0.25855 for the latter, as stated above. Each 
value (l/i')eff is then obtained from power-law fits to 
the FSS expression (|12(l . taking four consecutive points 
whose maximum size is L = Lmax- Results for the RDIS 
(empty circles) yield (l/i/)gg- = 1.478(5) at L^ax = 96. 
A rough estimate of the asymptote l/v, is obtained by 



extrapolating these points to the L~^^^ -^ axis, as 
seen in Fig. [T] (dotted hne). The extrapolation yields 
l/v « 1.464, well in agreement with previously reported 
results for the RDIS^^ Results for the AEIS (filled 
squares) give (l/j^)eff — 1.501(7) at L^ax — 96. Ef- 
fective exponents in this case clearly depart from values 
corresponding to the LRC fixed pointiiSi, through a region 
close but above the SRC fixed point at L^ax ~ 48. How- 
ever, at larger lattice sizes, greater values suggest that 
a another fixed point may rule the critical behaviour at 
the thermodynamic limit. 
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FIG. 7: (Color online) Effective correlation lenght inverse ex- 
ponent (l/f)gg- versus 1/Lmax, obtained by the FSS power- 
law fit l|12|l of four consecutive points ending at Lmax ■ Values 
for the RDIS (empty circles), extrapolated to Lmax -^ oo 
(dotted line) approach the result by Calabrese et al^. In 
the AEIS case (filled squares), points clearly depart from the 
LRC fixed point through a region close but above the SRC 
fixed point. Averages over disorder on less extensive simu- 
lations of the AEIS at the critical coupling Jc = 0.258570, 
give (l/i^)eff ^ 1-493 (diamond). The value (l/!^)eff ==; 1-497 
(triangle) is obtained from averaged reweighted curves at 
Jc = 0.258570. 

For the stable uncorrelated (SRC) disorder fixed point, 
the theorjii predicts that the exponent l/i^ should be 
smaller than 3/2. Additionally, the WH condition^: 
is well satisfied for this AEIS, where v is the pure 3DIS 
exponent, and a = 2{d — df) comes from LRC of gelling 
clusters (GC) within DLCA aerogels^. In effect, the frac- 
tal dimension for the GC within aerogels at c — 0.2 
is df « 2.2, as reported elsewhere^^. This condition, 
together with theoretical predictions reported by Prud- 
nikov et al}^, would give l/v « 1.4 for the 3DIS with 
LRC defects, at the corresponding a « 1.6. From Fig. [7| 
it is clear that (l/i^)eff is far above this value. Thus, it 
is not the LRC subset of disorder (the GC) which rules 
the critical behaviour of the AEIS, in the way it certainly 
does for 3DXY universality class in aerogels^. 

The exponent (l/z^)gfj at Lmax = 96 was also obtained 
taking averages of logarithmic derivatives from aver- 



age reweighted curves, at the critical coupling estimated 
above f Sect ion UTTEIl . J^ = 0.258570. For L = 56-96 we 
obtain (l/^)gg- w 1.497 (triangle, Fig. [Tj). In addition, 
less extensive additional realizations of the AEIS also for 
L = 56 — 96 (600 for each size), were made at this more 
accurate value Jc = 0.258570, and averages over disorder 
were taken directly from simulations. The power-law fit 
for these points gives an estimate (l/i/)gg- « 1.493 (di- 
amond). Although these effective exponents are lower 
than 3/2, there exists yet not enough evidence in this 
work that the RDIS fixed point would be reached at 



B. Specific heat and energy exponents 

To check our results about the correlation length ex- 
ponent for the AEIS, we study the FSS of the specific 
heat and the energy, at the simulation coupling Jc- 
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FIG. 8: ( To'p) Specific heat cn versus L for the AEIS in a 
linear-log scale. A logarithmic singularity (a = 0) for the 
specific heat at J^ follows from the fit (dashed line). This 
result reinforces the estimated Ijv « 1.5 (Fig|SJl. (Bottom) 
Energy plotted against i'""^)'", taking q = and 1/u — 3/2, 
and the corresponding linear fit is shown (dashed line). 

Top of Fig|Sl shows the specific heat Ch plotted versus 
L using a linear-log scale. The dashed line is a logarith- 
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mic fit using the four largest lattice sizes, 56 — 96. This 
results suggests that the singularity of the specific heat 
could be logarithmic, consistent with our result \/v « 1.5 
for L„iax = 96. After Josephson hyperscaling relation 
[a — 2 — dv) a specific heat exponent a « would be 
expected. Bottom of FigJHl shows the linear dependence 
of the energy E on _l("~i)/'^, taking v — 2/3 and a = 0, 
which confirms the results stated above. We made ad- 
ditional analisys to specific heat data, and the energy 
as well, using the scaling of both quantities in the case 
a/v < 0. This method was performed by Schultka and 
Manousakis in determining the (negative, very small) ex- 
ponent a/iy for the pure 3DXY model—. In the case 
a/i' < 0, the specific heat scales as Ch ~ Coo + cii"/'', 
while the energy scales as E = Eoo + -Eii^""^-*/". Non- 
linear fits to these expressions, using the six largest lattice 
sizes, 40 — 96, give self-consistent results a/i/ « —0.022 
and (a — l)/i^ ~ —1.512, in agreement with a correlation 
lenght exponent l/iy « 1.49. This result agrees with the 
tendency shown by the effective values in Fig[7| Exper- 
iments on the critical point of the LV transition of ''He 
in aerogelsii^ report a cusplike peak in the specific heat, 
but authors do not report an esimate for the exponent 
a. Our results point for the largest lattice sizes to a log- 
arithmic singularity, which may also be consistent with 
these experimental results. 
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FIG. 9: Order parameter (left) and susceptibility (right), at 
the simulated critical point Jc, versus L for the AEIS. The 
magnetization M (circles) scales as L~^''^ , and the averaged 
WolfF cluster size [(n.w)] scales aa— /,~2/3/i^ (squares), giving 
an average estimate of 2/3/i/ = 1.032(1). The FSS power-law 
fit for the susceptibility gives 'y/v — 2.044(4). All power- law 
fits have been made for L = 56 — 96. 



C. Magnetic exponents 

According to the FSS theory, the magnetization and 
the susceptibility scale as M ^ L^^/'^ and x ~ L'^^'^, 



respectively. In Fig. |21 we plot the order parameter M 
(left) and the susceptibility (right) as a function of L 
for the AEIS. Magnetization data (circles) have been fit 
to the preceding power-law FSS expression (continuous 
line), givieng /?/i/ — 0.523(3). Average sizes of Wolff 
clusters divided by L^, [('^w)] (squares), scale with the 
same exponent as the squared magnetization'*^. This is 
confirmed by the power-law fit (dashed line) which yields 
2/3/i^ == 1.019(6). These results give an average estimate 
2l3/iy = 1.032(6). Together with l/v w 1.5 this gives 
(3 = 0.34(4), close to the the pure 3DIS exponent and to 
the RDIS exponent. On the right side of Fig. O points 
for the susceptibility obtained from simulations near the 
critical point fit to the FSS power-law expression with 
the exponent 7/1/ = 2.044(4). All fits have been made 
for L = 56~ 96. 
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FIG. 10: (Color online) (a) Effective exponents for the sus- 
ceptibility give (7/i^)gfj- = 2 — r^gg- > 2. Extrapolation us- 
ing Lmax > 48 suggest that y/f < 2 in the thermodynamic 
limit, (b) Effective exponents (/3/i^)gg- for the the magne- 
tization (circles) and for the average mass of Wolff clusters 
(squares). At Lmax = 96 both approach the pure 3DIS the- 
oretical estimate— (f3/u)„g « 0.517. (c) The hyperscaling 
relation 7/1/ + 2/3/i^ = d — 3, not satisfied for effective expo- 
nents, tends to hold for larger lattice sizes. 

Results for magnetic effective FSS exponents, obtained 
by fitting four consecutive points from Fig. ending at 
Lmax, have been depicted in Fig. EH Fig. ITnT ai. shows 
the effective exponent for the susceptibility, (7/i^)effj 
plotted versus L^\^. This exponent tends to increase 
for Ljnax < 48, but beyond this size the tendency is to 
stabilize at a value close to that of percolation, "f/iy = 
2 — 77 « 2.05. For larger Lmax, it turns to approach a 
value lower than 2.00 (positive rj). The asymptotic ex- 



trapolated value seems to be the pure 3DIS exponent 
7/z/ = 1.966(3)^ or the RDISS^ 7/1/ = 1.963(5). Berche 
and collaborators^^ estimated effective values 7/t/ > 2 
for smaller concentrations in the RDIS case. Using the 
result for Lmax = 96, "f/v = 2.044(4), and our estimate 
V — 2/3, the exponent 7 found in the present work is 
1.363(9), slightly above the value 7 = 1.344(9) found in 
our MC simulations for the RDIS, and 7 = 1.342 ob- 
tained by Calabrese et alm^ for the same system. The 
tendency for larger L is to be closer to RDIS results. 

We must remind at this point, that theoretical most 
accurate results, by Prudnikov et al}^, predict a mag- 
netic exponent 77 < for the 3DIS in LRC disordered 
structures with an algebraic decay similar to that of the 
gelling clusters within DLCA aerogels at c = 0.2- . 

On Fig. ^1 (b) and (c) , squares represent effective ex- 
ponents (/3/i')eff obtained from [{n„)], and circles, those 
obtained directly from M. There is a strong variation 
of these results with Lmax- As stated above, the expo- 
nent obtained averaging both results using L„iax — 96 is 
l3/i' = 0.516(6), close to that of the pure 3DIS (3/i^ = 
0.517(3)2^. In addition, our results agree well with those 
reported for the RDIS case by Ballesteros et alm^, fi/v = 
0.519(3). In this work, we obtain P/v = 0.516(5) for the 
RDIS. The last effective value {(i/v)^^ [Lmax — 96), to- 
gether with 1/v « 1.5 gives an exponent /? = 0.343(9) 
for the order parameter. This result agrees well with ex- 
periments about the critical point of the LV transition of 
N2 in 95% porous aerogels (c = 0.05), reported by Wong 
et alM, which yield /3 = 0.35(5). 

Care must be taken with this agreement because, our 
results were obtained using DLCA at concentrations c — 
0.2, and these structures are quite different from those 
at c = 0.05. For the latter, most of impurities belong to 
the LRC gelling clusters, giving the DLCA aerogels a less 
random overall structure. 

As stated above, a possible explanation for these mag- 
netic exponents is the influence of the LRC disorder fixed 
point. The fractal dimension of the aerogel gelling clus- 
ter is df « 2.2, giving an exponent a ~ 1.6 associated to 
this structure^. Following Table IV from Prudnikov et 
alJ^, a value j/v > 2.0205 (77 < -0.0205) is expected. 
The effective value found in this work is close to this pre- 
diction, but it follows from Fig. Iior b) that a tendency 
exists to approach a value closer to the corresponding 
RDIS fixed point. Finally, using effective values the hy- 
perscaling relation •j/v + 2(3 /v = 3 seems not to hold, as 
seen in Fig. lior c). The violation of this hyperscaling rela- 
tion suggest that our results do not yet reach asymptotic 
values. Extensive simulations rest still to be performed 
at the more accurate value Jr = 0.258570. 



V. CONCLUDING REMARKS 

Extensive Monte Carlo simulations of the 3D Ising 
model (3DIS) with impurities have been reported in this 
paper. Using finite size scaling, critical couplings and 



exponents have been estimated for the 3DIS, in presence 
of randomly distributed impurities (RDIS), and confined 
in aerogel-like structures (AEIS). For the latter we have 
collocated Ising spins in the pores of simulated DLCA 
aerogels a.t c — 0.2. At this concentration, this objects 
are known to be non-fractal. However, the presence of 
hidden LRC could affect criticality, as predicted by the 
theoryiSi. It has been concluded elsewhersi that these 
LRC structures, the gelUing clusters, modify the critical 
behaviour of the 3DXY model, when confined in the same 
kind of aerogel- like structures. In the 3DIS case, however, 
our results for thermal exponents l/i^ ^ 1.5 and a/v ^ 
rest far above those for the LRC fixed point predicted by 
the theorjii^. Complementary simulations at a more ac- 
curate value of the critical coupling, Jc = 0.258570(13), 
give an exponent l/v ~ 1.49. Although similar thermal 
exponents have been reported by Pakhnin and Sokolo\i^ 
for the RDIS universality class, the asymptotic critical 
regime could have not been reached in our simulations, 
and more extensive simulations are yet to be performed 
at this more accurate Jc value. 

Effective critical exponents observed here for the AEIS 
change from a fixed point (LRC) at box sizes L < 48, to 
another (SRC) at box sizes L > 48 (Figs. rTlandlTU)!. 
probably indicating an oscillating approach to the stable 
fixed point. Theoretical predictions based on the Wein- 
rib and Halperin model^ii^, able to explain changes on 
the critical behaviour of the 3DXY model in the pores of 
DLCA aerogels^,, may also explain the influence of this 
type of disorder on the 3DIS. In this case, two competing 
effects are present: the random SRC subset of the disor- 
der (defined in Section Fill A|l . which already affects the 
critical behaviour of the 3DIS, and the LRC subset which, 
after the extended criterion (jSJ, may be relevant as well. 
This is certainly not the case for the 3DXY model, where 
Harris criterion prevents the SRC subset of impurities 
(islands) to be relevant: only the weak LRC distribution 
of impurities (GC) is relevant^ for the 3DXY. For the 
3DIS, theory predicts that both, LRC and SRC subsets 
may be relevant. Which one finally dominates the critical 
behaviour? 

Results presented in this paper suggest that, in the 
AEIS case, the critical behaviour is ruled by the SRC 
fixed point. A plausible explanation to these dominat- 
ing SRC effect is provided by theoretical work a^i^^ : for 
the 3DIS, RG flows converge to a more stable SRC fixed 
point, because at m = 1 the LRC fixed point is less 
stable (marginal). However, it has been mentioned be- 
fore, without a proof^, that amplitudes of disorder may 
in some cases affect criticality. Added in proof, we have 
to mention that in preliminary simulations of the 3DIS 
in presence of mixed kinds (LRC and SRC) of disorder—, 
evidence of a continuous flow from the LRC fixed point to 
the SRC one has been observed, when relative strengths 
are tunned from a pure LRC distribution of defects to 
a 1 : 1 proportion. In the AEIS case, we have analized 
relative amplitudes (strengths) of the LRC and the SRC 
subsets of disorder, for L = 128, to determine that up to 
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97% of defects are due to islands (SRC), while only 3% 
are due to the GC (LRC). 

To conclude, it has to be stated that the influence of 
aerogel-like distributions of impurities on the critical be- 
haviour of the 3DIS is yet far from being completely un- 
derstood. The problem is similar to that of phase tran- 
sitions in Ising systems with non-integer dimension or 
in fractal structures. It becomes clear that the fractal 
dimension, related to the exponent a of subjacent long- 
range correlations, may not be the only parameter to 
determine the universality class of the impure system^. 
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